Cross-correlation of CMB with large-scale structure: weak gravitational lensing 
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We present the results of a search for gravitational lensing of the cosmic microwave background 
(CMB) in cross-correlation with the projected density of luminous red galaxies (LRGs). The CMB 
lensing reconstruction is performed using the first year of Wilkinson Microwave Anisotropy Probe 
(WMAP) data, and the galaxy maps are obtained using the Sloan Digital Sky Survey (SDSS) 
imaging data. We find no detection of lensing; our constraint on the galaxy bias derived from the 
galaxy-convergence cross-spectrum is bg = 1.81 ± 1.92 {la, statistical), as compared to the expected 
result of bg ~ 1.8 for this sample. We discuss possible instrument-related systematic errors and 
show that the Galactic foregrounds are not important. We do not find any evidence for point source 
or thermal Sunyaev-Zel'dovich effect contamination. 

PACS numbers: 98.80.Es, 98.62.Sb, 98.62.Py 



I. INTRODUCTION 

Th e Wilkinson Microwave Anisotropy Probe (WMAP) 
|l08j | satellite has provided a wealth of information 
about the universe through its high-resolution, multi- 
frequency, all-sky maps of the cosmic microwave back- 
ground (CMB) Q. While the WMAP power spectrum 
[2 and temperature-polarization cross-spectrum Q are 
useful for probing the high-redshift universe (reioniza- 
tion and earlier epochs) 0, S S j the WMAP maps 
also provide an opportunity to study the low-redshift 
universe through secondary CMB anisotropics. While 
the effect of secondary anisotropics on the angular scales 
probed by WMAP {I < 700) is small compared to the 
primordial temperature fluctuations, the signal-to-noise 
ratio can be boosted by cross-correlating with tracers of 
the large scale structure (LSS) at low redshifts. Since 
the WMAP data release, several authors have used vari- 
ous tracers of LSS to measure the integrated Sachs- Wolfe 
(ISW) effect, the thermal Sunyaev-Zel'dovichftSZ) ef- 
fect, and microwave point sources la.l3. [liA[lU[l2.[l3 . [l4l| . 
The Sloan Digital Sky Survey (SDSS)Tloif is an excel- 
lent candidate for these cross-correlation studies due to 
the large solid angle covered at moderate depth. 

Another secondary anisotropy, which has not yet been 
investigated observationally, is weak lensing of the CMB 
by intervening large-scale structure. Weak lensing has 
attracted much attention recently as a means of directly 
measuring the matter power spectrum at low redshifts 
(e.g. _15|). The traditional approach is to use distant 
galaxies as the "sources" that are lensed to measure 
e.g . the matter power spectrum (e.g. jJa, [13j HM H^ 
12(1 |2]1 [2 ^. [23 1 ) or the galaxy-matter cross-correlation 
(e.g^lMISllM mill 113). However weak lensing 
of CMB offers an alternative method, free of intrinsic 
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alignments, uncertainties in the source redshift distribu- 
tion, and selection biases (since the CMB is a random 
field). Potential applications of CMB lensing described 
in the literature include precision measurement of cos- 
mological parameters ^3, 1^ 1^ ISS! and separation of 
the lensing contribution to the CMB B-mode polariza- 
tion from primoridal vector |34| and tensor perturba- 
tions [33, |3g, [33. While these applications are in the 
future, the WMAP data for the first time allows a search 
for weak lensing of the CMB in correlation with large- 
scale structure. This paper presents the results of such a 
search; our objective here is not precision cosmology, but 
rather to detect and characterize any systematic effects 
that contaminate the lensing signal at the level of the 
current data. This step is a prerequisite to future inves- 
tigations that will demand tighter control of systematics. 

In this paper, we perform cross-correlation analy- 
sis between the CMB weak lensing field derived from 
WMAP and a photometrically selected sample of lu- 
minous red galaxies (LRGs) in the SDSS at redshifts 
0.2 < z < 0.7. The photometric LRGs are well-suited for 
cross-correlation studies because of their high intrinsic 
luminosity (compared to normal galaxies), which allows 
them to be observed at large distances; their high num- 
ber density, which suppresses shot noise in the maps; and 
their uniform colors which allow for accurate photometric 
redshifts and hence determination of the redshift distri- 
bution. We use the measured cross-spectrum between the 
lensing field and the projected galaxy density to estimate 
the LRG bias bg. At the present stage, we are using the 
bias as a proxy for the strength of the cross-correlation 
signal, just as has been done in recent analyses of the 
ISW effect |l[l3,[ll[ll[l3|; we are not yet trying to use 
the bias in cosmological parameter estimation, although 
this is a possible future application of the methodology. 
We do not have a detection of a cross-correlation, and 
hence our measured bias bg = 1.81 ± 1.92 is consistent 
with zero. 

This paper is organized as follows. The most inrpor- 



tant aspects (for this analysis) of the WMAP and SDSS 
data sets, and the construction of the LRG catalog, are 
described in Sec.m The theory of CMB lensing and re- 
construction methodology are explained in Sec. IIIII The 
cross-correlation methodology and simulations are cov- 
ered in Sec. IIVI and the results are presented in Sec. 
We investigate possible systematic errors in Sec. lVIl and 
conclude in Sec. I VI II Appendix ^ describes the spher- 
ical harmonic transform algorithms and associated con- 
ventions used in this paper, and Appendix \n\ describes 
the algorithm used for the C~ operations that arise in 
our analysis. 



issues, namely the loss of frequency-dependent informa- 
tion (useful as a test of foregrounds), the inability to 
separate cross-correlations between different DAs from 
auto-correlations (useful to avoid the need for noise bias 
subtraction), and the complicated inter-pixel noise cor- 
relations (due to the smoothing used to create the map 
and the varying weights of the different frequencies). 
The foreground-cleaned map of Ref. |43| recovers the full 
WMAP resolution, but the practical difficulties (for the 
purpose of lensing reconstruction) associated with ILC 
still apply. We have not used either of these maps in this 
paper. 



II. DATA 



B. LRG density from SDSS 



A. CMB temperature from WMAP 

The WMAP mission ^38] is designed to produce all-sky 
maps of the CMB at multipoles up to I ^several hundred. 
This analysis uses the first public release of WMAP data, 
consisting of one year of observations from the Sun-Earth 
L2 Lagrange point. WMAP carries ten differencing as- 
semblies (DAs), each of which measures the difference in 
intensity of the CMB at two points on the sky; a CMB 
map is buit up from these temperature differences as the 
satellite rotates. ( WMAP has polarization sensitivity but 
this is not used in the present analysis.) The DAs are 
designated Kl, Kal, Ql, Q2, VI, V2, Wl, W2, W3, and 
W4; the letters indicate the frequency band to which a 
particular DA is sensitive [l|, y^ (the K, Ka, Q, V, and 
W bands correspond to central frequencies of 23, 33, 41, 
61, and 94 GHz, respectively). The WMAP team has 
pixelized the data from each DA in the HEALPix [llOJ 
pixelization system at resolution 9 0, |43 . This system 
has 3,145,728 pixels, each of solid angle 47.2 sq. arcmin. 
These maps are not beam-deconvolved; this, combined 
with the WMAP scan strategy, results in nearly uncorre- 
lated Gaussian uncertainties on the temperature in each 
pixel. 

In this paper, we only use the three high-frequency 
microwave bands (Q, V, and W) because the K and Ka 
bands are very heavily contaminated by galactic fore- 
grounds and have poor resolution. (The foreground 
emission is not a Gaussian field and cannot be reliably 
simulated, so in cases where it dominates over CMB 
anisotropy and instrument noise, we cannot compute reli- 
able error bars on the cross-correlation.) For the galaxy- 
lensing correlation, we have used the sky maps produced 
by the eight high-frequency DAs. The variances of the 
temperature measurements are obtained from the effec- 
tive number of observations Nobs ■ 

Note that the WMAP "internal linear combination" 
(ILC) map [l| cannot be used for lensing studies because 
of its degraded resolution (1 degree full width half maxi- 
mum, FWHM), which eliminates the multipoles I ~ 350 
of greatest importance for the lensing analysis. An ILC- 
based lensing analysis would also suffer from practical 



The Sloan Digital Sky Survey (SDSS) [43 is an ongoing 
survey to image approximately n steradians of the sky, 
and follow up approxima tely one million of the detected 
objects spectroscopically |43,I3. The imaging is carried 
out by drift-scanning the sky in photometric conditions 
|45| . in five bands (uariz) 149 . l47l | using a specially de- 
signed wide-field camera [43 . These imaging data are 
the source of the LSS sample that we use in this paper. 
In addition, objects are targeted for spectroscopy using 
these data 49] and are observed with a 640-fiber spec- 
trograph on the same telescope. All of these data are 
processed by completely automated pipelines that detect 
and measure photometric properties of objects, and as- 
trometrically calibrate the data [50, Isj]. The SDSS is 
well under way , and has had three major data releases 
[53,|53i|53,|53; this paper uses all data observed through 
Fall 2003 (296,872 HEALPix resolution 9 pixels, or 3,893 
square degrees). 

The SDSS detects many extragalactic objects that 
could, in principle, be used for cross-correlation with sec- 
ondary anisotropics ^3- The usefulness of LRGs as a 
cosmological probe has been appreciated by a number of 
authors j57i, i5^ . These are typically the most luminous 
galaxies in the universe, and therefore probe cosmologi- 
cally interesting volumes. In addition, these galaxies are 
generically old stellar systems and have extremely uni- 
form spectral energy distributions (SEDs), characterized 
only by a strong discontinuity at 4000 A. The combi- 
nation of these two characteristics make them an ideal 
candidate for photometric redshift algorithms, with red- 
shift accuracies of Uz ^ 0.03 |53. We briefly outline the 
construction of the photometric LRG sample used in this 
paper below, and defer a detailed discussion of the selec- 
tion criteria and properties of the sample to a later paper 

Our selection criteria are derived from those described 
in Ref. 58J. However, since we are working with a photo- 
metric sample, we are able to relax the apparent luminos- 
ity constraints imposed there to ensure good throughput 
on the SDSS spectrographs. We select LRGs by choos- 
ing galaxies that both have colors consistent with an old 
stellar population, as well as absolute luminosities greater 



than a chosen threshold. The first criterion is simple to 
implement since the uniform SEDs of LRGs imply that 
they lie on an extremely tight locus in the space of galaxy 
colors; we simply select all galaxies that lie close to that 
locus. More specifically, we can define three (not inde- 
pendent) colors that describe this locus, 



c_L = (r - i) - 0.25(5 -r)- 0.18 , 
d± = {r - i) - 0.125{g - r) , 
CM EE 0.7(g-r) + 1.2(r-i-0.18) 



(1) 



where g, r, and i are the SDSS model magnitudes ,52j in 
the g, r and i bands (centered at 469, 617, and 748 nm re- 
spectively). We now make the following color selections. 



Cut I: 

Cut II: d. 



|ci|<0.2; 
> 0.55, g - r > 1.4. 




(2) 



Making two cuts (Cut I and Cut II) is convenient since 
the LRG color locus changes direction sharply as the 4000 
A break redshifts from the g to the r band; this division 
divides the sample into a low redshift (Cut I, z < 0.4) 
and high redshift (Cut II, z > 0.4) sample. 

In order to implement the absolute magnitude cut, 
we follow 1^ and impose a cut in the galaxy color- 
magnitude space. The specific cuts we use are 



Cut I : rpetro < 13.6 + — , rpetro < 19.7, 
Cut II: i< 18.3 + 2rfi, i < 20, 



(3) 



where rpetro is the SDSS r band Petrosian magnitude 
[53 • Finally, we reject all objects that resemble the point- 
spread function of the telescope, or if they have colors 
inconsistent with normal galaxies; these cuts attempt to 
remove interloping stars. 

Applying these selection criteria to the ~ 5500 degress 
of photometric SDSS imaging in the Galactic North 
yields a catalog of approximately 900,000 galaxies. Ap- 
plying the single template fitting photometric redshift al- 
gorithm of [53, we restrict this catalog to galaxies with 
0.2 < Zphoto < 0.6, leaving us with ~ 650,000 galaxies. 
We use the regularized inversion method of Ref. 59] as 
well as the photometric redshift error distribution pre- 
sented there, to estimate the true redshift distribution 
of the sample. The results, comparing the photometric 
and true redshift distributions are shown in Fig. ^ Fi- 
nally, this catalog is pixelized as a number overdensity, 
g — 5n/n^ onto a HEALPix pixelization of the sphere, 
with 3,145,728 pixels. We also mask regions around 
stars from the Tycho astrometric catalog 61], as the 
photometric catalogs are incomplete near bright stars. 
The final catalog covers a solid angle of 3,893 square de- 
grees (296,872 HEALPix resolution 9 pixels) and contains 
503,944 galaxies at a mean density of 1.70 galaxies per 
pixel. 



FIG. 1: The LRG redshift distribution. The black histogram 
shows the photo-2 distribution, the red curve is the true red- 
shift distribution estimated by regularized deconvolution of 
the photo-z errors. 



III. LENSING OF CMB 



A. Definitions 



Gravitational lensing re-maps the primordial CMB 
anisotropy T into a lensed temperature T according to 



r(n) = f(n + d(n)). 



(4) 



where the 2-vector d is the deflection angle of null 
geodesies. To first order in the metric perturbations, d 
can be expressed as the gradient of a scalar lensing po- 
tential, d = V$, where V is the derivative on the unit 
(celestial) sphere. We may also define the convergence 
K = — iy • d. Assuming the primordial CMB is statis- 
tically isotropic with some power spectrum Ci , it can be 
shown J63, |63 that the multipolc moments of the lensed 
temperature field have covariancc 
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(5) 



where we have introduced the Wigner 3j symbol, and the 
coupling coefficient is 



Juih 



{2L + l)i2h + l){2l2 + 1) 



L{L + l)\l 167r 

X {[L(i + 1) + h{h + 1) - h(h + l)]C/i 

+ [L{L + 1) - h{h + 1) + hih + 1)]^ } 

h h L 




(6) 



The WMAP satellite does not directly measure T, but 
rather a beam-convolved temperature: 



r"(n) = / B"(n,n')T(n')d2n' + noise. 



(7) 



(Here a=Kl..W4 are the differencing assemblies.) In 
most of this analysis we have approximated the beam 
by the WMAP circularized beam transfer function JGI]. 
The temperature multipole moments recovered assuming 
a circular beam are 



rpOL 






(8) 



where Bf are the beam transfer functions. If the beam 
is truly circular, Eq. IjSJl returns an unbiased estimator of 
the beam-deconvolved CMB temperature; in Sec. IVICI 
we will consider the effect of the WMAP beam ellipticity 
on lensing estimation. Note further that T^ is only well- 
determined up to some maximum multipole / because the 
B" drop to zero at high I. 



B. Theoretical predictions for lensing 



In this paper, we aim to measure the galaxy- 
convergence cross-correlation, Cf , where g = 5n/n 



is the projected fractional overdensity of galaxies; this 
section briefly presents the theoretical prediction from 
the ACDM cosmology. In a spatially flat Friedmann- 
Robertson- Walker universe described by general relativ- 
ity, the convergence is given in terms of the fractional 
density perturbation 5 by: 



K(n) = '^ttGnPq 



XiXCMB - x) 
XCMB 



{l + z)6ix,n)dx, (9) 



where x is the comoving radial distance, z is the redshift 
observed at radial distance x, Po is the present-day mean 
density of the universe, and xcmb is the comoving dis- 
tance to the CMB. The galaxy overdensity does not come 
from a "clean" theoretical prediction, but on large scales 
it can be approximated by 



5(ft) = 



JbgAf{x)S{x,-n)dx 
S^{x)dx ' 



(10) 



where the matter power spectrum P^ is evaluated at co- 
moving wavenumber k = l/x and at the redshift corre- 
sponding to conformal time 770 — X- It is obtained us- 
ing the transfer functions from CMBFast (35| and the 
best-fit 6-parameter fiat ACDM cosmological model from 
WMAP and SDSS data from Ref. [H (f^^/i^ = 0.0232; 
n^K^ = 0.1454; h = 0.695; r = 0.124; erg = 0.917; 
TT-s = 0.977). We have found that varying each of these 
parameters over their \a uncertainty ranges gives a ±18% 
effect for erg (for which Cf scales as oc cr|) and ± < 14% 
effect for the other parameters. Since the overall sig- 
nificance of the detection is only 0.9cr, this dependence 
of the template Cf on cosmological parameters will be 
neglected here. 

The lensing signal is on large scales and so we have not 
used a nonlinear mapping. The peak of the LRG redshift 
distribution is at z ^ 0.5, corresponding to a comoving 
angular diameter distance ^ \.'Ah~^ Gpc, in which case 
the smallest angular scales we use [l = 300) correspond to 
k = 0.23ft. Mpc"i and A^[k) = 0.7. The nonhnear evolu- 
tion at this scale according to the Peacock-Dodds formula 
[63 is a 10% correction to the matter power spectrum 
and is thus much smaller than the error bars presented 
in this paper (although it is not obvious what this im- 
plies about the galaxy-matter cross-spectrum, which is 
the quantity that should appear in Eq. 1111) . Future ap- 
plications of CMB lensing in precision cosmology will of 
course require accurate treatment of the nonlinear evo- 
lution. 

In this paper, we will assume that galaxy bias bg is 
constant so that it may be pulled out of the integrals in 
Eqs. H10() and (|ll|l . If the bias varies with redshift (as 
suggested by e.g. ^68]), then the best-fit value of bg will 
be some weighted average of bg over the redshift distri- 
bution; this need not be the same weighted average that 
one observes from the auto-power spectrum, since the 
latter is weighted differently. Computation of the auto- 
power in photometric redshift slices |60| suggests that 
over the redshift range 0.2 < z < 0.6 the bias varies from 
1.7-1.9; this variation can safely be neglected given our 
current statistical errors. Note however that a detection 
of 6g 7^ via a fit to Eq. (|ll|l assuming constant bg would 
rule out Cf'*' = and hence would be robust evidence for 
a galaxy-convergence correlation, regardless of the red- 
shift dependence of the bias. 



where Af{x) is the distribution in comoving distance and 
bg is the galaxy bias. (The SDSS LRG sample is at low 
redshift z < 0.7 and does not have a steep luminosity 
function at the faint end of our sample, so we neglect the 
magnification bias.) For I ^ 1 and smooth power spectra 
for matter and galaxies, this may be approximated by the 
Limber integral: 
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, (11) 



C. Lensing reconstruction 

We construct lensing deflection ma ps u sin g qu adratic 
reconstruction methods [ISiliilzaEilzlllElzlIzi, 

which have been shown to be near-optimal for lensing 
studies of the CMB temperature on large (/ < 3500) 
scales 176]. Non-quadratic methods may be superior if 
CMB polarization is used [73, or on very small scales 
Izl El US; these cases are not of interest for WMAP, 
since the sensitivity is insufficient to map the CMB polar- 
ization and arcminute scales are unresolved. Quadratic 



estimation takes advantage of the cross-coupling of differ- 
ent multipoles induced by gravitational lensing, namely 
the 0{klm) term in Eq. (jHJ. The maximum signal-to- 
noisc statistic for CMB weak lensing is the divergence 
of the temperature- weighted gradient. We construct, 
for each pair of differencing assemblies a and f3, the 
temperature- weighted gradient vector field G"^: 



G°'l^{n) = U[wf"]{h)W[CWf'^]{n) 



-[Wf^]{n)V[CWf"]{n) 



(12) 



where WT and CWt are defined by the convolution re- 
lations [Wf%^ = Wif^^ and [CWf%^ = Cmfl'^. 
Note that Eq. p2|l is exactly the same as the G statis- 
tic of Ref. j73| except that we have multiple differenc- 
ing assemblies, and we have left open the choice of the 
weight Wi. While Wi = {d + C,"°*«e)-i jg statistically 
optimal, there are also practical considerations that af- 
fect this choice. Specifically, it is desirable that the WT 
and CWT convolutions are almost-local functions of the 
CMB temperature (to minimize leakage from the Galac- 
tic Plane), and that the same Wi be used for all differenc- 
ing assemblies (so that any frequency dependence of our 
results can be attributed to foregrounds or noise, rather 
than merely a change in which primary CMB modes we 
are studying). We choose the following weight function: 



Wi 



Ci + (0.0346lAiK^)e 



2^„K'+l)/300^ 



(13) 



which clearly has the optimal Cf dependence in the high 
signal-to-noise regime. Note that in the range of I we use 
{I < 800), the Wi drop to zero with increasing I faster 
than the Q-band beam transfer functions. Hence the 
computation of WT'^^'^^ are stable even though T'^^'Q^ 
are beam-deconvolved. Due to their narrower beams, 
this stability also applies to the V and W-band DAs. 
The Kl and Kal DAs have wider beams and hence lens- 
ing reconstruction using the weight Eq. (|13|) is unsta- 
ble for these DAs. We set Wo = Wi = to reject the 
monopole (not observed by WMAP) and dipole signals. 
The power spectrum Ci used for the lensing reconstruc- 
tion is WMAP best-fit ACDM model with scalar spectral 
runninga^ to the CMB data ( WMAP-^ ACBAR 
+CBI |S2|)- Errors in the Ci used in the lensing recon- 
struction cannot produce a spurious galaxy-temperature 
correlation because they result only in a calibration error 
in the lensing estimator. Furthermore, WMAP has de- 
termined the Ci to within several percent (except at the 
low multipoles, which give a subdominant contribution 
to both WT and WCWT), whereas the lensing cross- 
correlation signal is only present at the la level, this 
error is not important for the present analysis. 

In the reconstruction method of Ref. [73, a filtered 
divergence of G is taken to extract the lensing field. We 
avoid this step because it is highly nonlocal, and hence 



can smear Galactic plane contamination into the regions 
of sky used for lensing analysis. In principle, we would 
prefer to directly cross-correlate G with the LRG map, 
but this too is difficult because the noise power spectrum 
of G is extremely blue. Wc compromise by computing v, 
a Gaussian- filtered version of G: 



„(II.-L) 



-l{l + l)afj/2M\\, 



(14) 



Here || and _L represent the longitudinal (vector) and 
transverse (axial) multipoles, which are the vector ana- 
logues of the tensor E and B multipoles. The Gaussian 
filter eliminates the troublesome high-Z power present in 
G(n) and makes v suitable for cross-correlation studies. 
We have chosen a width ctq = 0.01 radians (34 arcmin). 
The vector field v can be written in terms of the tem- 
peratures T directly in harmonic space. The longitudinal 
components are given by 
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(15) 



where we have defined: 
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We will not need the formula for the transverse compo- 
nents. While we perform the reconstruction (Eg s. 1 121 and 
I14|l in real space, the harmonic-space relation fEa. ll5|l is 
useful for computing the response of the estimator, and 
for estimating foreground contamination and beam ef- 
fects. In particular, from the orthonormality relations 
for Wigner 3j-symbols, we have 



„"/3(ll)\ 
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21 ^\ 



l^lm = RlK,lr. 



(17) 



which defines the calibration of v as an estimator of 
the lensing field. The response factor Ri is shown in 
Fig. 12 we have verified this response factor in simula- 
tions fSec UvTl . 

There are 36 pairs of differencing assemblies that could 
be used to produce estimated lensing maps v^^^: the 8 
"autocorrelations" (a = (3) and 28 "cross-correlations" 
(a ^ P). Note that instrument noise that is non- 
uniform across the sky can produce a bias in the 
"autocorrelation" -derived lensing map v"". In principle 
the noise bias could be estimated and subtracted, just as 
can be done for the power spectrum. However, this is 
dangerous if the noise properties are not very well mod- 
eled. Since the WMAP noise is in fact strongly variable 
across the sky we only use the "cross-correlation" a ^ (3 
maps. 



Lensing response factor 



(a) Lensing weights in harmonic space 
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FIG. 2 
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Multipole, L 
The response factor Ri of Eq. 1171 satisfying 



(«L"*"^) - Rl-^r. 



One problem that we find with this method is that 
the V field contains "ghosts" caused by the Galactic 
plane (where small-scale temperature fluctuations of sev- 
eral millikelvin or more can occur due to Galactic emis- 
sion). We solve this problem by setting T = within 
the WMAP Kp4 jg^j Galactic plane mask. We have veri- 
fied that using the Kp2 mask instead produces only small 
changes to the results. 

The weight functions Wi and CiWi are shown in Fig.O 
We also show the real-space weights, given by 



2^ + 1, 



W^(0) = ^^— VK/P^cos 



(18) 
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(b) Lensing weights in real space 
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FIG. 3: (a) The weight functions Wi and CiWi. (b) The same 
weight functions in real-space (Eg. I18II : WT and CWT are 
obtained by convolving the (beam-deconvolved) temperature 
T with these kernels. 



and similarly for [CbF](0). 



D. Frequency-averaged lensing maps 

The methodology outlined in Sec. If If Gl allows us to 
construct 28 lensing maps v"^ corresponding to the 28 
pairs of differencing assemblies. For this analysis, we 
need to produce an "averaged" lensing map v^-^-^^ based 
on a minimum- variance linear combination of the 28 DA- 
pair maps. The averaged lensing map is determined by 

ITT) 

the weights a^o : 



Q/3 



(19) 



We select these weights to minimize the amount of power 
in v(^^\ subject to the restriction X)^a/3 1^ this is 

done by minimizing the total vector power in v between 

multipoles 50 and 125: V = Ei=5o EL=-, k;^™"^^ "P, 

(TT) 

which is a quadratic function of the weights a^ ^ . The 
optimal weights a^^ are complicated to establish ana- 
lytically since the maps v"'^ are highly correlated. We 
have therefore minimized V using a simulated lensing 
map (see Sec lIVC)) . Using a simidated map rather than 





-' —xSr-^^ 



FIG. 4: The divergence of the lensing vector field map, 
V • v' , smoothed with a 30 arcmin FWHM Gaussian, and 
displayed in Galactic MoUeweide projection. Note the promi- 
nent artifacts surrounding the Galactic plane cut and the 
point sources (which are removed by the Kp05nS10\ps2 cut). 



the real data avoids the undesirable possibility of the 
weights being statistically correlated with the data. We 
also fix Oq^ Q2~^ because the ■y'^^'Q'^ map would be the 
most heavily contaminated by point sources. The weights 
so obtained are shown in Tabled A map of V • v*^"'"-'"-', 
smoothed to 30 arcmin resolution (Gaussian FWHM) is 
shown in Fig. 01 

Also, to study foreground effects on lensing estima- 



TABLE I: The weights Ua/B used for the 
lensing map a^^ , and the six "individual 



overaU-averaged 
frequency" maps 



JQQ) 



(WW) 



DA pair (a/3) 


ATT) 
"a/3 


„{"^l''2) 


Q1,Q2 


0.000000 


1.000000 


Q1,V1 


0.063759 


0.220605 


Q1,V2 


0.075370 


0.267859 


Q2,V1 


0.061538 


0.227341 


Q2,V2 


0.098884 


0.284196 


Q1,W1 


0.031583 


0.167730 


Q1,W2 


0.026896 


0.113943 


Q1,W3 


0.018391 


0.094086 


Q1,W4 


0.028571 


0.130429 


Q2,W1 


0.023816 


0.158925 


Q2,W2 


0.017153 


0.107379 


Q2,W3 


0.014494 


0.097062 


Q2,W4 


0.027164 


0.130444 


V1,V2 


0.106330 


1.000000 


V1,W1 


0.053278 


0.150880 


V1,W2 


0.029381 


0.094494 


V1,W3 


0.029664 


0.087425 


V1,W4 


0.035463 


0.107193 


V2,W1 


0.048598 


0.168035 


V2,W2 


0.049979 


0.137365 


V2,W3 


0.034621 


0.115438 


V2,W4 


0.046029 


0.139169 


W1,W2 


0.016453 


0.199003 


W1,W3 


0.014448 


0.170408 


W1,W4 


0.019908 


0.232958 


W2,W3 


0.014514 


0.129381 


W2,W4 


0.003365 


0.133767 


W3,W4 


0.010347 


0.134484 



tion, we would like to construct averaged lensing maps 
y(QQ) ^ T^iQ"^) ^ etc. where we only average over differenc- 
ing assemblies at the same frequency, thereby preserving 
frequency-dependent information. There are six of these 
maps (QQ, QV, QW, VV, VW, and WW); the last col- 
umn of Table^shows the weights used to construct them. 
The power spectrum of the longitudinal mode of v ob- 
tained on the cut sky (Kp05nS10\pS2 cut, which excludes 
point sources; see Sec. lIVAp is shown in Fig. |S1 



Template projection is dangerous for 
cross-correlation studies involving galaxies because the 
dust template of Ref. [SJI is used to extinction-correct the 
LRG magnitudes, thus template errors could introduce 
spurious correlations between the CMB and galaxy maps. 
Since visual inspection of the uncleaned WMAP maps re- 
veals Galactic contamination outside the Kp2-rejected re- 
gion at all five frequencies, we have used the more conser- 
vative KpO mask in Sec. I VI Dl for our galaxy-temperature 
correlations. Because the SDSS covers only a small por- 
tion of the sky, we speed up the cross-correlation compu- 
tation by using only WMAP data in the vicinity of the 
SDSS survey region. We define the "SIO" region to con- 
sist of those pixels within 10 degrees of the SDSS survey 
area. The KpOnSlO cut accepts 774,534 HEALPix pixels 
(10,157 sq. deg.). 

When analyzing primary CMB anisotropics, it is cus- 
tomary to mask detected point sources in order to elim- 
inate this spurious contribution to the temperature. For 
secondary anisotropy studies, the analysis should be 
done both with and without the point sources because 
the point sources may correlate with large scale struc- 
ture, hence naively masking them could lead to mis- 
leading results. Therefore for the galaxy-temperature 
correlation used in Sec. IVIDI we have constructed the 
KpOnS10\ps mask by rejecting all pixels within the 
WMAP point source mask (with a 0.6 degree exclusion 
radius around each source). The KpOnS10\ps cut accepts 
756,078 HEALPix pixels (9915 sq. deg.). 

For the lensing analysis, we must use a more conserva- 
tive mask than KpO because the lensing estimator v is a 
nonlocal function of the CMB temperature, hence v(n) 
responds to foreground emission several degrees away 
from n. We have therefore constructed a "Kp05" mask 
consisting of all pixels within KpO that are at least 5 de- 
grees away from the KpO boundary; the Kp05nS10 mask 
used for the lensing analysis accepts 753,242 HEALPix 
pixels (9,878 sq. deg.). We have also constructed a point 
source- removed version, Kp05nS10\pS2 , in which all pix- 
els within 2 degrees of the point sources are rejected. This 
mask accepts 598,795 HEALPix pixels (7,853 sq. deg.). 



IV. CROSS-CORRELATION COMPUTATION 

A. Sky cuts 

In some regions of the sky, particularly the Galactic 
plane, microwave emission from within the Milky Way 
and from nearby galaxies dominates over the cosmolog- 
ical signal. For their CMB analysis, the WMAP team 
removed this signal by (i) masking out a region based on 
a smoothed contour of the K-band temperature, which 
they denote "Kp2" [83], and (ii) projecting out of their 
map microwave emission templates for synchrotron, free- 
free, and dust emission based on other observations 



B. Galaxy-convergence correlation 

Having constructed the vector field v, we proceed to 
compute its cross-correlation Cf^ with the LRG map. 
We construct the data vector 



i^LRG^^lens) 



(20) 



jiLRG) 



of length iVv^ -I- 2Ny^ ' , where y^LRG is a vector 



jiCMB) 
pix ' 

containing the galaxy overdensities g — 6n/n in each 
SDSS pixel, and x/e„s consists of the two components of 
V at each WMAP pixel. (We will suppress the frequency 
indices {QQ), (QV), etc. on v for clarity; it is under- 
stood that the analysis below is repeated for each pair of 
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FIG. 5: The longitudinal mode power spectrum C"jj within the Kp05nS10\ps2 cut (solid line), for each of the six frequency 
pairs. The points show 10 simulations (containing no lensing or foregrounds). Since the purpose of this plot is to compare the 
simulations to the actual data, the sky cut has not been deconvolved, rather we have plotted the average of | J v-Y;^ d^n//afeyp 
within bands of width Al = 10. 



frequencies.) The covariance of x is then: 



C EE (XX^) = 



Q{LRG) Q{x) 



q{x)T Qilens) I ■ (21) 

The cross-correlation matrix C^^^ has components: 

^ITk = E ^r>l™(nOY|L(n,) • e^, (22) 



Im 



where i represents an SDSS pixel index, j is a 
WMAP pixel index, and K = 9,(1) indicates which com- 
ponent of the vector v is under consideration. We bin 
the cross-spectrum Cf" into bands. 






A / 1 ImhM) < I < ^max(^) 



otherwise 



(23) 



and take the c as the parameters to be estimated. 

In order to construct an optimal estimator for the 
galaxy-convergence cross-spectrum, we need a prior auto- 
correlation matrix for the LRGs and for the lensing map. 
(This is a "prior" in the sense of quadratic estimation the- 
ory J83, ISfl I93, Isll , and has nothing to do with Baycsian 
priors.) We take a prior of the form 



C 



(LRG) 



Im 



CfYC^{n,)YUi^, 



NS^. 



(24) 



where Cf ® is the galaxy power spectrum (excluding Pois- 
son noise) and N is the noise variance per pixel. We 
have taken A'' to be the reciprocal of the mean number 
of galaxies per pixel, appropriate for Poisson noise (we 
use the mean galaxy density per pixel to avoid biases 
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FIG. 6: The prior power spectrum used for the LRG map, 
obtained by apphcation of a pseudo-C; estimator to the SDSS 
scan region. The dashed hue shows the Poisson noise. 



associated with preferential weighting of pixels contain- 
ing fewer galaxies). The prior power spectrum C; 
is determined by application of a pseudo-C; estimator to 
the LRG maps; the resulting power spectrum is shown 
in Fig. El We have set Cf^ = 0.01 > C^l^ for Z = 0, 1 
to reject the galaxy monopole ("integral constraint") and 
dipole modes from the cross-correlation analysis. 

It can be shown [SS, 90, 91] that for Gaussian data 
with small C^^-*, the optimal estimator for the c^ would 
be obtained by taking the unbiased linear combinations 
of the x|;^cC('^^'^)-iPAC('^"'')-ix,ens. Thcsc estima- 
tors are frequently called "QML" or quadratic maximunr- 
likelihood estimators, although they are, strictly speak- 
ing, not maximum-likelihood. We do not have the full 
matrix C'''^"''', and our only knowledge of this matrix 
comes from simulations. Thus we have instead con- 
structed, for each lensing map v"'', the quadratic com- 
binations: 






'-LRG^ 



dc^ 



X;e 



(25) 



where A represents a band index. This differs from the 
QML estimators in that the C^^ weighting is applied 
only to the LRGs, while uniform weighting is applied to 
the CMB lensing map v; thus Eq. H25|) can be viewed 
as a sort of half-QML, half-pseudo-C; estimator for the 
cross-spectrum. The expectation value of Qa can be de- 
termined from Eq. (|21|l : it is: 



(Q. 



Tr 



dcB 



q(LRG) 



_idC^' 



dc 



c^ = Rabc^, 



(26) 
which defines the response matrix Rab ■ Note that unlike 
the response matrix of the optimal quadratic estimator, 
Rab is not equal to the Fisher matrix. The trace in 
Eq. if^ may be computed using a stochastic-trace algo- 
rithm: 



i?^^ = ((C(^«G)-iy): 



x)T 



,ac(^)9c( 



dc^ dc^ 



(27) 



where y is a random vector of length n!^^^ ' consist- 
ing of ±1 entries. The vectors in Eqs. (P5|l and (|?7jl are 
constructed in pixel space. Harmonic space is used only 
as an intermediate step in the convolutions required to 

compute the matrix-vector multiplications, e.g. q^b y; 
these are computed by the usual method of converting to 
harmonic space, multiplying by dCf^/dc^, and convert- 
ing back to real space. The matrix inverse operations are 
performed iteratively as described in Appendix IbI This 
method allows us to easily compute estimators for the 
band cross-powers. 



[R 



'ItAB 



Qi 



(28) 



While this estimator is manifestly unbiased, we do not 
know its uncertainty because we do not know the co- 
variance C''*^"*^ of the v-field. We determine the uncer- 
tainty via a Monte Carlo method: we construct random 
CMB realizations according to the null hypothesis of no 
lensing in all eight DAs used for the lensing reconstruc- 
tion, and feed them through the lensing reconstruction 
pipeline fSec. lIIIC|) . 

In order to estimate the galaxy bias from the binned 



cross-power spectrum estimators c 
the response of each estimator c^ 
d{c^)/dbg. This is given by 



, we need to know 
to the galaxy bias, 



d{c 



dbn 



= [R-'] 



-UAB 



Tr 



dbn 



q{LRG) 



_iac( 



X) 



dc^ 



(29) 



which is computed by a stochastic-trace algorithm anal- 
ogous to Eq. H27|l . The galaxy-v correlation matrix 
dC(^)^/d6g 



IS 



.(X) 



'^^CjK ^y^R^'^^'^ 



dbn 



Im 



dbn 



Yr^{n,)Yl{h,) ■ BK 



(30) 



where Ri is the lensing response factor of Eq. (|17|l . If we 
knew C^^ , it would be optimal to use C^^ ^^ weighting, 
in which case we could simply use dC^^^^ /dbg as a cross- 
power template with no loss of information. However 
since we have not calculated C^^, and our only informa- 
tion on this covariance matrix comes from the ability to 
generate random realizations of v, we cannot do this. 



C. Simulations 

Simulating lensed and unlensed CMB maps is neces- 
sary both for verifying the analysis pipeline as well as 
for determining the optimal weighting of Sec. IIII Dl The 
general procedure used here is to generate a simulated 
primary CMB temperature Tim, convergence Kim, and 
galaxy density fluctuation Snim in harmonic space. These 
are Gaussian random fields and hence it is a straightfor- 
ward matter to produce random realizations from the 
power spectra and cross-spectra of T, k, and Sn. After 
generation of the realization, the primary temperature 



10 



and deflection field (generated from the convergence, as- 
suming an irrotational deflection field) are pixelized in 
HEALPix resolution 10 (12,582,912 pixels of solid an- 
gle 11.8 arcmin^ each). The lensed CMB temperature T 
is then computed in real space from Eq. I^J. Because 
the "deflected" HEALPix pixels no longer lie on curves 
of constant latitude, we use a non-isolatitude spherical 
harmonic transform (see Appendix fXl we have used pa- 
rameters L' ~ 6144 and K — \1 since high accuracy is 
required) to evaluate Eq. Q. The beam convolution 
relevant to each DA is then applied by converting to har- 
monic space, multiplying by Bf and the pixel window 
function, and converting back to real space. Finally the 
simulated CMB temperature field T'(n) is degraded to 
HEALPix resolution 9, and appropriate Gaussian "in- 
strument" noise is added independently to each pixel. 
Note that the resolution 10 pixels are used here to im- 
prove the fidelity of the simulation, in particular to ensure 
that the effects of the elongated HEALPix pixels on the 
lensing estimator are properly simulated. 

A crude model for the WMAP beam ellipticity is in- 
corporated into the simulations as follows. At each point, 
we have 



^"(")=9 E Y.Bi,TimYC,{n){e 



ifj.^p\ 



(31) 



side— A,_B Im^ 



where the Y^^ are spin-weighted spherical harmonics and 
the beam moments B;^ arc the multipolc moments of the 
beam in instrument-fixed coordinates (with the "North 
Pole" along the boresight and the (p' — meridian in the 
scan direction): 



Bif^ — 



An 
21 + 1 



B{n,„st)Yi*{n,nst)d'^n,, 



(32) 



The average value {e^^"^) is taken over the position an- 
gles of the instrument when n is scanned. The sum over 
sides is over the two sides of WMAP. Equation H31|) is 
only an approximation because (i) the two sides of the 
differencing assembly may not scan each pixel exactly 
the same number of times, or may have slightly differ- 
ent weights; and (ii) because WMAP is a differential in- 
strument, {T" {9 , (j))) is also affected by beam-ellipticity 
effects on other parts of the sky. Since the two beams of 
a given DA are separated by ~ 140 degrees, this results 
in an "echo" of a g iven microwave source at a separation 
of 140 degrees [4fl (and higher-order echoes should also 
be present); we have neglected these. 

We have made several further approximations to 
Eq. H31() in order to speed up the simulations. First, we 
have only included the ellipticity modes fi = ±2, since 
these dominate the difference between the azimuthally 
symmetrized beam and the true beam. The beam ellip- 
ticity is thus described by the real and imaginary parts 
of Bi^2', recall B/,_2 = B*2- We have calculated Bi^2 by 
taking the non-isolatitude spherical harmonic transform 
of the WMAP beam maps ^j. Secondly, because the 



side A and B beams are approximate mirror-images of 
each other, we have only considered the component of 
the beam ellipticity along the scan direction. The com- 
ponent of the beam ellipticity at 45 degrees to the scan 
direction is suppressed because, to the extent that the 
side A and B beams are mirror images and scan each 
pixel the same number of times, this component cancels 
in Eg. 13 II when we sum over the two sides. 

Finally, we have used a simple model for the scan pat- 
tern (e*'"'^) for each DA. The WMAP scan pattern is 
crudely approximated as a rotation around the space- 
craft —Z axis, followed by a precession of this axis in a 
22.5 degree radius circle around the anti-Sun point, fol- 
lowed by rotation of the anti-Sun point along the ecliptic 
plane. Relative to the spacecraft —Z axis, the effective 
number of observations in one rotation is: Nobs{0,4>) — 



-Nobs{< 



.2iip\ — 



KS{coa6 — cos6'c)/27r, where dc is the an- 



gle between the instrument boresight and the spacecraft 
— Z axis, and K is a. constant. We can convert these to 
harmonic space in the spacecraft coordinates, 

[Nob.]im{-Z) = KSrnoYioiecO); 
[Nobs{e''^)]i^{^Z) = -K5rnoYiUec,0). (33) 

Averaged over the precession cycle of WMAP, this be- 
comes 



[A^ofcs];o(av) 
[iVo,,(e2^^)],„(av) 



i^FKcos22.5 deg)Pi{0)Yio{9c,0y, 
-KPi{cos22.5 deg) 
xPi{O)Yil{0,,O). (34) 



(The 771 7^ moments vanish.) Once these have been 
obtained, we may transform back to real space and find 
{e^^"^) by division. This is then rotated from the ecliptic 
to the Galactic coordinate system. To speed up com- 
putation, the elliptical correction to the beam was only 
computed on the HEALPix resolution 9 grid whereas the 
dominant circular part was computed on the resolution 
10 grid and then degraded by pixel-averaging. 



V. RESULTS 
A. Galaxy-convergence cross-spectrum 

The individual cross-spectra obtained at different fre- 
quencies are shown in Fig. [T] The frequency-averaged 
cross-spectrum is shown in Fig. |S1 both with and with- 
out point sources. 



B. Amplitude determination 



We estimate the bias amplitude hg by fitting the ob- 
served galaxy-convergence cross-spectrum C[''^ to the 
theoretical model, Eq. 111|) . We begin by obtaining the 
covariance matrix F of c as determined from M = 50 
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FIG. 7: The galaxy-convergence correlation using the Kp05nS10\ps2 cut, for each of the six combinations of frequencies. The 
error bars are strongly correlated across different frequencies. The dashed curve shows the theoretical signal for our best-fit 



value of bg = 1.81. 



simulations: 



\AB 



M 



M 



(^ 



'm 



(35) 



(The simulated c^ are generated by producing a random 
realization of the CMB as described in Sec. IIVCI with 
no lensing, feeding it through the lensing pipehne, and 
correlating it against the real SDSS LRG map.) The bias 
is estimated from the weighted average of the observed 



cross-powers in our A^ = 14 bands: 



bo = 



2Za=i c 



A d<^th if-AA 
db„ /-^ 



1^A= 



i($) /r^^' 



(36) 



where c^ is the theoretical prediction for the binned 
(j23|l : this is directly proportional to 
The response dcf^/dbn is obtained 



cross-spectrum, Eq. 
the bias, 



^th 



from Eq. H29I) . It is trivially seen that Eq. (|36|l is an 
unbiased estimator of bg, regardless of the covariance of 
the c . Since we are working in harmonic space with 
bands of width AZ > 20 > A9~^, where A0 is the typi- 
cal width of the survey region (in radians), the different 
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FIG. 8: The galaxy-convergence correlation, using the 
frequency- averaged v'^"^' map. The top panel shows the cor- 
relation using the Kp05nS10\ps2 mask (which rejects point 
sources). The middle panel shows the correlation using the 
Kp05nS10 mask, which does not reject point sources. The 
bottom panel shows the correlation with the v'"^"^^ field ro- 
tated by 90 degrees; this should be zero in the absence of 
systematics (see Sec. I VI ^ . The error bars are from simula- 
tions as described in Sec. IV Bl The dashed curve shows the 
theoretical signal for bg — 1.81. 



Z-bands are very weakly correlated, so we have not at- 
tempted to further optimize the relative weights of the 
various cross-power estimators c^. 

The most obvious way to estimate the uncertainty in bg 
by noting that Eq. \6<5\\ is a hnear function of the c^ , and 
substituting in the covariance matrix T^^ of the {c^}: 



uibg] incorrect) 



E 



N 
A,B=1 



Y^B dcth dct^ /pAApSB 
dhr, dhfj ' 



ELt ($) /r- 



(37) 

This calculation is incorrect for finite number M < cxd 
of simulations because it neglects the fact that the T^^ 
are themselves random variables. One approach to the 
problem is to take a sufficiently large number of simu- 
lations M that the error in Eq. H37|l becomes negligible. 
The difficulties in this approach are that it could be very 
computationally intensive; we do not know whether M 
simulations are "enough" unless we try even larger val- 
ues of M to check convergence. An alternative method, 
which we have used here, is to run an additional M' = 50 
simulated realizations of {c^} (identical to those used to 
compute F"^^ except for the random number generator 
seed) , compute the bias hg from them, and then compute 
their sample variance. The resulting error bars can be 
analyzed using the well-known Student's t-distribution. 
The "Icr" t error bars (which have 49 degrees of free- 
dom) obtained by this method are shown in Table ^ 
The mean bias values obtained from these 50 random re- 
alizations are shown in the "random" column in the table. 
Also shown in Table ITTl (in the "foreground" columns) are 
the results obtained by feeding the Galactic foreground 
templates of Sec. IVI E"2l through the lensing pipeline and 
correlating these with the real LRG map. 

In the case of the Kp05nS10 cut (last column in Ta- 
ble]^, which does not reject point sources, the v^'^'^', 
v('3^)^ and v^'^M') maps have power spectra that are 
boosted significantly by point source contamination (see 
Sec. IVID l|l . Therefore, even if the correlation of the 
point sources with the galaxies can be neglected, the er- 
ror a{bg) obtained in these bands for the Kp05nS10 mask 
is probably underestimated, as noted in the table. 

The x^ values for fits to zero signal are shown in Ta- 
ble lllll These are obtained using the 14 band cross-power 
spectra (Fig.|SJl, and the 14 x 14 covariance matrix is ob- 
tained from 100 simulations, 



X' - [f-i (100 Sims)] Ai3C 



^c^. 



(38) 



Because the number of simulations is finite, there remains 
some noise in this covariance matrix and this must be 
taken into account in interpreting the x^- In particular, 
the x^ variable does not exactly follow the standard x^ 
(so we have denoted it with a hat). The distribution 
and p- values can, however be calculated as described in 
Ref. [93, Appendix D. As noted previously, the errors 
for the Kp05nS10\pS2 mask in the QQ, QV, and QW 
combinations are suspect. 
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TABLE II: The bias bg estimated using several frequency ranges (first column). The second column gives the bias bg obtained 
with the point sources removed (Kp05nSfO\ps2 cut); these are the numbers that should be thought of as our "result." The 
uncertainties are la with t errors (49 d.o.f.); the first error given used the elliptical beam simulations, the error in parentheses is 
obtained from circular beam simulations (Sec. lVICl . The "TT" frequency combination is a weighted average of QV, QW, VV, 
VW, and WW. The column labeled "foreground" show the bias bg obtained by correlating the v maps derived from the Galactic 
foreground maps of Sec. I Vf E2l against the LRG map. Similarly, the column labeled "random" show the bias bg obtained using 
the average of the 50 random realizations of v in place of the WMAP-derived v map. The column labeled "TT weight" shows 
the bias determined by using the same weights as a function of I as for the TT frequency combination in Eq. 1361 : of course 
this has no effect on bg{TT). The final column, labeled "KpOSnSlO ", is the bias obtained without the point source cut. 



Frequency 


Bias, bg 


foreground 


random 


TT weight 


bg, Kp05nS10 


QQ 


+3.62±4.48(±4.33) 


-0.001640 


-0.34 ±0.63 


+3.13 + 4.55 


+6.30 + 3.93" 


QV 


+3.10±2.19(±2.00) 


-0.000001 


+0.07 ±0.31 


+3.39 + 2.29 


+3.93 + 1.99" 


QW 


-0.11±2.53(±2.35) 


-0.000480 


+0.06 ±0.36 


-0.20 + 2.56 


+0.48 + 2.43" 


VV 


-1.41±2.95(±3.13) 


-0.000737 


+0.88 + 0.42 


+0.34 + 3.04 


+0.58 + 2.35 


VW 


+2.63±2.56(±2.11) 


-0.000617 


+0.35 + 0.36 


+2.58 ± 2.62 


+2.51 + 2.03 


WW 


+0.23±3.11(±2.71) 


-0.000817 


+0.20 + 0.44 


+0.65 ± 3.22 


-0.75 + 2.63 


TT 


+1.81±1.92(±1.72) 


-0.000449 


+0.23 + 0.27 


+1.81 + 1.92 


+2.43 + 1.58" 



"Because of point sources, the v maps from these bands contain 
excess power in the Kp05nS10 region. Thus the simulation error 
bars shown here are likely underestimates. 



TABLE III: The x^ values obtained for fits to zero signal 
from the galaxy-convergence cross-spectrum. The first col- 
umn shows the results with the Kp05nS10 mask; the sec- 
ond using the Kp05nS10\ps2 mask; and the third using the 
Kp05nS10\ps2 mask with 90 degree rotation of v. The y^ has 
14 degrees of freedom (the 14 Z-bins shown in Fig. |SJ and the 
covariance matrices were obtained from the 100 simulations 
described in Sec. IV Bl As described in the text, the finite 
number of simulations means that the expectation value of 
the x^ is not 14 but is larger due to uncertainty in the covari- 
ance matrix; the mean of the x^ distribution is 16.5 and the 
standard deviation is 6.8. We have also given the cumulative 
probability distributions. 



Freq. 



Kp05nS10 



Kp05nS10\ps2 

X 



■'' P{< x') 



Kp05nS10\ps2 
+90 degrees 

Pi<x'] 



X 



QQ 


51.51" 


0.9996" 


16.24 


0.55 


28.24 


0.940 


QV 


37.01" 


0.990" 


27.51 


0.931 


21.89 


0.81 


QW 


16.74" 


0.58" 


11.74 


0.26 


21.58 


0.80 


VV 


19.08 


0.70 


11.86 


0.27 


24.85 


0.89 


VW 


25.30 


0.90 


20.28 


0.75 


20.06 


0.74 


WW 


11.54 


0.25 


11.36 


0.24 


21.82 


0.81 


TT 


26.90" 


0.922" 


19.97 


0.74 


30.71 


0.963 



"Because of point sources, the v maps from these bands contain 
excess power in the Kp05nS10 region. Thus the simulation error 
covariance matrices are likely underestimates, and hence the x^ 
and P{< x^) values are suspect. 



VI. SYSTEMATIC ERRORS 



E and B modes, and in the absence of systematics there 
should be no B-mode signal. In the case of CMB lens- 
ing using the vector estimator v, the analogous test is to 
rotate v by 90 degrees (thereby interchanging the longi- 
tudinal and transverse parity modes). This rotated map 
can be fed through the cross-correlation pipehne in place 
of the original v. In the absence of systematics, this gives 
zero signal; the error bars need not be the same as for the 
longitudinal modes, but they can still be determined from 
simulations as described in Sec. lVBl The cross-spectrum 
is shown in Fig. |S1 and the x^ values in Table IIIII 

The lowest-Z point in the rotated cross-spectrum 
(Fig. |SJ is 3.4cr negative. It is difficult to assess the sig- 
nificance of this anomaly since it is an a posteriori de- 
tection {p — 0.00127 for the two-tailed i-distribution) ; in 
any case, it is responsible for the relatively high x^ value 
(p = 0.037) in the Kp05nS10\pS2 +90 degree column of 

Table UTTl It is unlikely that this correlation {g*v\^ ) 
represents any real astrophysical or cosmological effect, 
since it violates parity. This anomaly is also distinct 
from the much-discussed "low quadrupole" observed by 
WMAP, since the former is based on a high-pass filtered 
CMB map with power coming predominately from CMB 
modes with I ^fewxlO^. Another possible explanation 
would be some source of excess power in the v map at 
low I, which would increase the error bar relative to sim- 
ulations and thus lower the statistical significance of this 
point. However if we take the v maps, and compute the 
un-deconvolved power spectrum 



A. Ninety-degree rotation test 

One of the standard systematics tests in weak lensing 
studies using galaxies as sources has been to rotate all 
of the galaxies by 45 degrees and look for a shear signal. 
The 45 degree rotation is used because it interconverts 



p=i:i: 



/<20 m=-l 



E 

neKp05nS10\ps, 



v(fi).y,;^;(n) 



(39) 



for both the real v map and the 100 simulated maps, we 
find that the real map has the 28th highest value of P out 
of 101 maps, i.e. there is no evidence for excess power. If 
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this point is due to some systematic, it must be present 
at all three frequencies, since this point is negative by at 
least Icr in all of the frequency combinations except QQ, 
where the binned Cf^ from the rotated map at ^ < 20 is 
(1.0 ±3.5) X 10-^ 

It is thus difficult to explain the lowest-/ point in Fig.|S| 
point in terms of any systematic error. The true test for 
whether this is in fact just a statistical fluctuation is to 
wait for the error bars to become smaller with future 
WMAP data and see whether this point becomes more 
significant or goes away, and in the former case whether 
it exhibits a frequency dependence. 



TABLE IV: The calibration factors obtained via end-to-end 
simulations, for the Kp05nS10\ps2 mask and various combi- 
nations of frequencies. Error bars are la, t-distributed with 
49 degrees of freedom. 



Frequencies 


Calibration factor. 


1±C 


QQ 


1.17±0.13 




QV 


1.24 ±0.09 




QW 


1.19 ±0.07 




VV 


1.31 ±0.11 




VW 


1.15 ±0.07 




WW 


0.92 ±0.10 




TT 


1.18 ±0.07 





B. End-to-end simulation 

Another important systematic test is to verify, in an 
"end-to-end" simulation, that the lensing estimator and 
Cf^ cross-spectrum estimator are calibrated properly. 
This can be done as follows. We run 50 simulations 
in which simulated Gaussian g and k maps are gener- 
ated with the cross-spectrum appropriate for bg = 1. 
The K map has the power spectrum C;"" expected for a 
ACDM cosmology, while the g map is constructed from 
gim = Ci'^Kim/Cj^^. In principle one could add addi- 
tional noise to g to boost its power spectrum to match 
the observed Cf^, but there is no reason to do this as 
it increases the number of simulations required and has 
no effect on the calibration. The k maps are then used 
to generate lensed CMB maps using the simulation code 
described in Sec. lIVC^ and the output temperature maps 
Tqi...Tw4: fed through the lensing reconstruction pipeline 
and then the Cf^ estimator. Finally, we estimate the bias 
in each simulation using Eq. I|36() . This output bg is the 
calibration factor 1 ± C appropriate for cross-correlation 
studies. 

The calibration factors obtained from this procedure 
are shown in Table IIVI the table reveals that the lens- 
ing pipeline is calibrated at the ^ 20% level. Calibra- 
tion factors of this order have been observed in previous 
simulations ITo. l94l and have been investigated analyti- 
cally [7^ 123, ISa, ISgl , where the main effect has been the 
non-linear lensing effects (i.e. the order k^ and higher 
terms that have been dropped in the Taylor expansion, 
Eq. Isj. In our case, the calibration error may also have 
a contribution from the elliptical beam. In any case, the 
calibration errors are not significant at the level of the 
current data (i.e. no detection). 



C. Beam effects 

We have used a crude model for the WMAP beam 
ellipticity. An incorrect model for the beam can have 
three effects on the galaxy-convergence cross spectrum 
and hence on the bias determination: it can (i) produce 
a shear calibration bias in the v map (which may depend 
on the wavenumber I and orientation of the convergence 



mode in question, and may vary across the sky as the 
effective beam varies); (ii) modify the noise covariance 
matrix of v; and (iii) introduce artifacts (i.e. biases) in 
the V map because it invalidates the assumption that the 
signal is statistically isotropic. The calibration problem 
is obviously of concern for attempts to do precision cos- 
mology with lensing. However since we do not have a 
detection, the only effect of the calibration bias is to af- 
fect our upper limits on the lensing signal. The change 
in the noise covariance of the lensing map v is poten- 
tially more serious because it can alter the variance of our 
cross-spectrum estimator and hence affect the statistical 
significance of any lensing detection. Because artifacts in 
the V map do not correlate with the galaxy distribution, 
they are essentially also a source of spurious power, and 
for cross-correlation measurements they are only a con- 
cern if their power spectrum is comparable to that of the 
noise. 

The effect of the beam on the noise covariance can be 
addressed as follows. We have re-computed the uncer- 
tainties a{bg) (see Sec. IVBjl using 50 simulations with 
a circular beam instead of our elliptical beam model. 
The uncertainties are shown in parentheses in Table ^1 
They are at most 20% different from the error esti- 
mates obtained from the elliptical beams (but this may 
not be significant because the cr(6g) values from simula- 
tions are themselves drawn from a random distribution 
- namely, the square root of a x^ distribution with 49 
degrees of freedom - and hence have an uncertainty of 
l/-\/2 • 49 « 10%). Since replacing our model of the beam 
ellipticity with the inferior model of a circular beam has 
only a < 20% effect on (j{bg), it is doubtful that (j{bg) 
would be altered by more than this by use of an improved 
beam model. 

Finally, we come to the issue of the calibration. Our 
estimator for the bias was constructed assuming a circu- 
lar beam, and given that the true beam is not circular, we 
expect that it may be mis-calibrated, i.e. (bg) = {1+C)bg, 
where C is the calibration bias. At present, the best 
way to test for such a bias is via simulations, such as 
those of Sec. IVIBI There we found a calibration bias of 
C = 0.18 ± 0.07, which is not important at the level of 
the present data. 
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D. Extragalactic foregrounds 

The lensing estimator of Eq. (|14() will respond not only 
to real lensing signals, but to any other perturbations 
of the CMB. Of greatest concern is the contamination 
from extragalactic foregrounds, which may induce spuri- 
ous correlation of the lensing estimator with the galaxy 
distribution since the extragalactic foregrounds (tSZ and 
point sources) are expected to correlate with large scale 
structure. The presence of the extragalactic foregrounds 
causes the observed temperature T(n) to be incremented 
by some amount AT"(n). Assuming that the extra- 
galactic foregrounds are not correlated with the primary 
CMB or with instrument noise, this causes the expecta- 
tion value of V (averaged over primary CMB and noise 
realizations) to be incremented by (compare to Eq. I15|l 
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Thus the possible source of contamination of the 
convergence-galaxy correlation signal Cf^ is the corre- 
lation of the LRGs with quadratic combinations of the 
foreground temperature. In the case of the tSZ fore- 
ground, the contribution to Eq. (|40|l can be broken up 
into a "single-halo" term in which the two factors of 
AT come from the same halo, and a "two-halo" term 
in which the two factors of AT come from different ha- 
los. The single-halo term exists even if the tSZ halos are 
Poisson-distributed, whereas the two-halo term acquires 
a nonzero value only from clustering of the halos. Much 
of this section will be devoted to an investigation of the 
properties of the quadratic combinations in Eq. (|^ and 
an assessment of their magnitude. Unfortunately, we will 
see that this does not result in useful constraints on the 
foreground contamination to our measurement of hg. 



1. Point sources 

It is readily seen that point sources are a major 
contribution to the power spectrum of v, especially 
in the lower- frequency bands. This can be seen from 
Fig. 1^1 in which the power spectrum C'^^, is shown in 

the Kp05nS10 region (in which point sources are not 
masked). The v^'^'^^ v^'^^)^ and v^'^^) maps are heavily 
contaminated while for the higher-frequency maps point 
sources are subdominant to CMB fluctuations. While the 
contribution to the v autopower in Kp05nS10 is large, we 
are interested here in whether - and at what frequencies 
- the point source contribution to v"'^ correlates with the 
LRG map when the Kp05nS10\pS2 mask (which masks 
point sources) is used. 

A single point source with frequency-dependent flux 
Fa{i') (in units of blackbody /iK sr) at position Aa will 



produce a spurious contribution to the temperature of 

ATr^^Fa{v^)Y:^{ha). (41) 

Plugging this into Eq. ()4U|I . we find that the shift in (v) 



is: 



A(«r™^"^> = Fa(a)^a(/3)r,;.(n.)rps(0, (42) 



where the response function rps{l) is: 



rpsil) 



' {21' + 1){21" + 1) 
4tt{21 + 1) 



I I' I" 




/C 



ii'i" 



(43) 



In Table El we show the product Fa{i'i)Fa{i'2) for sev- 
eral possible point source spectra. Note that for the steep 
spectra characteristic of WMAP point sources (a ~ 0.0), 
the contamination of v^'^Q) and v^Q^^ is far greater than 
contamination of the higher-frequency lensing maps. 
Therefore these lower-frequency bands are a useful test 
of point source contamination of the galaxy-convergence 
correlation. The dependence of the estimated bg on the 
combination of frequencies, Fa{h'i)Fa{i'2), is shown in 
Fig. 1101 If there were point source contamination of our 
measurement with spectral index a = 0.0, the points 
in Fig. ^1 would be expected to fall roughly along a line 
bg ex Fa (h'l) Fa {1^2)', this is only rough because the weight- 
ing of different I bins is slightly different at different fre- 
quencies in Eq. H36|l . and the contaminating signal need 
not have the same angular dependence as the galaxy- 
convergence correlation. If we re-calculate the six fre- 
quency combinations bg{QQ)...bg{WW) using the same 
weighting of different I bins as for the bgiTT) measure- 
ment, we get the values in the column in Table Hll labeled 
"TT weight." Assuming a synchrotron-like spectrum for 
the point sources, a correlated least-squares fit of the 
form 



bg{viV2) = 6f 



S ^ ' {FaFa){TT) 



(44) 



to the various frequency combinations will return for 
bg a point-source-marginalized measurement of the bias, 
and for Abg {TT) a measurement of the point-source 
contamination to the unmarginalized bg{TT). The re- 



l(0) 



[PS), 



suits of such a fit are h^g' = 0.58±2.36 and Ahf' {TT) = 
0.73 ± 1.18; the two measurements are of course anti- 
correlated with correlation coefficient p ~ —0.63. There 
is thus no evidence for point source contamination, al- 
though the statistical errors are too large to definitively 
say whether such contamination is present at the level of 
the signal. It would be useful to have lower-frequency in- 
formation here in order to improve the constraints, how- 
ever this is not possible as the CMB multipoles used in 
our analysis {I ^ 350) are not resolved by WMAP K- and 
Ka-band differencing assemblies. 

The contamination in Cf^ from a point source can also 
be estimated from angular information since the point 
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FIG. 9: The cut-sky power spectrum of the longitudinal mode of v, C,"m; the power spectrum of the actual map is shown 
as the solid line, while the points are from 10 simulations. This figure is identical to Fig. |3 except that we have used the 
Kp05nS10 cut (i.e. in this figure the point sources are not masked), and we have not shown the VW or WW spectra (which do 
not differ significantly from Kp05nS10\ps ). Because of contamination by point sources, the power in the QQ and QV maps 
is greater than in the simulations or the point source-cut maps. 



sources have an angular dependence (roughly Poisson) 
unlike that of the CMB. Equations (glj and (g2l give 



AC; 



gK{TT) 
I 



rpsjl) 
Ri 



y^iTT) F{a)F{p) 



'■aP 



afi 



F{Q? 



xF{Q)ACt 



'gTiQ) 



(45) 



where AC; 



'ffT(Q) 



is the point source-induced galaxy- 



temperature cross spectrum in Q-band. If take a typi- 
cal spectrum a = 0.0 and a flux F{Q) = 1 Jy (i.e. the 
brightest point sources not excluded by the WMAP point 
source mask |83,|93), we find 



V- iTT) F{a)F{(3) 
F{QY 



ap 



F{Q) = 6.7 X IO^VK sr. (46) 



The ratio rps{l)/Ri is plotted in Fig. ^] Of course, 
not all of the point sources have F{Q) = 1 Jy, but this 
is the worst-case scenario since ACf^ oc FACf , hence 
if the galaxy-temperature cross-spectrum is coming from 
fainter sources the contamination will be even less. (This 
scaling with F occurs because the spurious contribution 
to the galaxy-convergence cross spectrum is quadratic 
in the flux whereas the contribution to the galaxy- 
temperature cross spectrum is linear.) We have com- 



puted the Q-band galax y-te mperature cross-spectrum us- 
ing a QML estimator [Qg on the KpOflSlO cut; if we 
take this cross-spectrum, and assume that at Z > 60 the 
cross-spectrum is entirely due to point sources with the 
synchrotron spectrum, the derived contamination to the 
galaxy-convergence spectrum is as shown in Fig. I12r a). 
One can also use the difference between Q and W-band 
cross-spectra; in this case, if a synchrotron spectrum with 
a — 0.0 for the point sources is assumed, the difference 
Cf^iQ) - Cf^{W) must be multiphed by 1.295 to re- 
cover the Cf (Q) used in Eq. H45|l . This result is shown 
in Fig. I12r b) ; the error bars (obtained from simulations) 
are now smaller at low I because the CMB fluctuations 
are suppressed. The point-source induced error in the 
bias can be computed from the data in Fig. IT^ b) by 
plugging these ACf^ values into Eq. (|SS|l : the result is 

Abf'^^ = -0.14 ±0.51. 



We conclude that the point-source contamination to 
Cf is at most of the same order as the signal in this range 
of multipoles. If one ignores correlations between distinct 
point sources so that Eq. H45|l is valid and assumes the 
a = 0.0 spectrum, then Fig. ll2r b) suggests that the point 
source contamination is less than the observed signal. 
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Frequency dependence of computed galaxy bias 
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FIG. 10: The dependence of the computed galaxy bias bg on 
frequency. The horizontal axis is Fa{i'i)Fa{i^2) /0.2'i9, com- 
puted for a typical synchrotron spectrum (q — 0.0); the nor- 
malization is chosen so that the frequency-averaged TT result 
has Fa (i^i) Fa {1^2)/ 0.249 = 1. The horizontal line is zero, and 
the point with dashed error bars is TT. (This point is of 
course the average of the other data points and contains no 
additional information.) The error bars are correlated; the 
X^ for a frequency-independent bias is 6.31 with 5 degrees of 
freedom. 



Lensing estimator response to point sources 
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FIG. 11: The ratio rps/Ri describing the response of the lens- 
ing estimator v to contamination from Poisson-distributed 
point sources. Note that on large scales this is negative, i.e. 
regions with more point sources are interpreted by the lens- 
ing estimator as regions of negative convergence (underdense 
regions). At high I the situation is reversed. 



2. Thermal Sunyaev-Zel'dovich effect 

In the case of the tSZ effect, the frequency depen- 
dence is exactly known; see Table Unfortunately 
this frequency dependence is extremely weak in the 
WMAP bands, with Fa{i'i)Fa{i'2) varying by only a fac- 
tor of 0.667 from the QQ to WW bands. Therefore we 
must resort to the angular dependence to separate tSZ 
from lensing. 

In order to perform an analysis similar to that of 
Sec. IVID II we need to determine or place a bound on 
the galaxy-SZ power cross-spectrum, and determine the 
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FIG. 12: (a) The contamination to the galaxy-convergence 
cross-spectrum due to point sources, based on the Q-band 
galaxy-temperature cross-spectrum in the KpQnSlO region 
(see text). The dashed line is the best-fit signal from Sec. lV B1 
with bg = 1.81. We have removed the I < 60 points since these 
may contain ISW signal, (b) The same figure, except derived 
from the difference of galaxy-temperature cross spectra in Q 
and W bands. This cancels any contribution from ISW and 
reduces cosmic- variance errors at low I. 



maximum flux \F\ of a tSZ source (we use |_F| since tSZ 
sources have negative net flux in the WMAP bands). 
In principle, tSZ haloes can be extended, however since 
our lensing estimator uses information from multipoles 
I ^ 350 (physical wavenumber k/a = 0.7 Mpc"^ at 
Dap > 0.5 Gpc), the haloes will not be resolved. This 
argument of course only applies to the single-halo contri- 
bution to Eq. H40|l . The flux from a tSZ source is 



F = -fbT, 



CMB 



CTT ksiTe, 



M 



n2 



= -8 X 10"'*/iK sr 



10i6M(T, keV 



X 

X coth — 
2 



0.18 



(DapX 



4 — a; coth — 
2 



(47) 



where ut is the Thomson cross secton, fc^ is Boltzmann's 



constant, /{, = ^b/^r. 



and 



are the electron and 



proton masses, (Te) is the mean electron temperature, 
fie'Tip is the baryonic mass per free electron. Dap is the 
physical angular diameter distance, M is the total mass 
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TABLE V: The point source frequency-dependent response 
Fa{i^i)Fa{v2) for several point source spectra. Included are 
three power-law spectra, Tantenna oc v"'^ (a — —1.0, 0.0, and 
-1-0.5) and a nonrelativistic tSZ spectrum. We have normal- 
ized to unit response in v''^'^'. The "TT" frequency combina- 
tion is an average over the other pairs of frequencies weighted 
by a^J"^^p. 



Bands 




Fa{Vl)F, 


.i'^2) 






Q = -1.0 


a = 0.0 


a = 4-0.5 


tsz 


QQ 


1.000 


1.000 


1.000 


1.000 


QV 


0.320 


0.476 


0.581 


0.946 


QW 


0.099 


0.228 


0.345 


0.817 


VV 


0.102 


0.227 


0.338 


0.895 


VW 


0.032 


0.109 


0.200 


0.773 


WW 


0.010 


0.052 


0.119 


0.667 


TT 


0.137 


0.249 


0.350 


0.838 



of the halo, and x — hv/ksTcMB = v jhl GHz. The 
frequency-dependent factor 4 — x coth(x/2) ranges from 
1.91 (Q band) to 1.57 (W band). For tSZ sources that are 
physically associated with the LRGs (distance z > 0.2 or 
Dap > 0.5 Gpc) we will have |i^| < 6 x 10~^fj,K sr even 
for extremely massive {M{Te) = IQ-'^^M© keV) haloes. 
If we estimate contamination to the galaxy-convergence 
correlation in analogy to Eq. (|45|l . we find 
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(48) 



if we take |F| = 6 x 10"'^ /iK sr, the limits on the contam- 
ination from the Q-band correlation Cf are similar to 
the limits for point sources from Cf (Q) (cf. Fig. 112b.). 
Unfortunately, like our similar analysis for point sources, 
this analysis of the tSZ contamination does not tell us 
anything new, since if the contamination from tSZ in our 
bg measurement were large compared to the statistical 
error of ±1.92, we would have measured the wrong bg. 
Worse, it only applies to the single-halo contribution in 
Eq. H4UII . whereas the two-halo contribution is likely to 
dominate on sufficiently large scales. 



E. Galactic foregrounds 

Foreground microwave emission from our own Galaxy 
can introduce spurious features in the weak lensing map. 
Because of their Galactic origin, these features cannot 
be correlated with the LRG distribution. However, it is 
possible that they can correlate with systematic errors 
in the LRG maps, most notably (i) stellar contamina- 
tion of the LRG catalog and (ii) incomplete correction 
(or over-correction) for dust extinction. We have used 
two methods to address these potential problems. The 
first (Sec. IVI E l\i is to correlate the CMB lensing map 
V with stellar density and reddening maps. The sec- 
ond (Sec. IVI E 2p is to correlate the LRG density map 
with simulated lensing contamination maps Av obtained 



TABLE VI: The bias obtained by substituting stellar den- 
sity {5ni,/ni,) and extinction [E{B — V)] maps in place of the 
LRG map (Sng/fig), for the Kp05nS10\ps2 sky cut. The fi- 
nal column, A6g, is the error in the galaxy bias if we use the 
coefficients c* and ce from the linear fit (Eg. 1491 to estimate 
the contamination of the LRG map by stars and extinction. 



Frequency 


stars 


E{B - V) 


Abg 


QQ 


-14.38 


-0.17 


+0.02 


QV 


-h2.37 


+0.19 


+0.10 


QW 


-2.12 


+0.01 


+0.03 


VV 


+3.40 


+0.13 


+0.05 


VW 


-2.29 


+0.14 


+0.11 


WW 


+3.99 


+0.24 


+0.11 


TT 


+0.71 


+0.15 


+0.09 



by feeding microwave foreground templates through the 
lensing pipeline. 



1. Stellar density and reddening tests 

The dominant systematics in the SDSS that could cor- 
relate with Galactic microwave foregrounds are stellar 
contamination of the LRG catalog and dust extinction. 
We study these by constructing two maps: a map from 
SDSS of the density (5n^/n^ of "stars" (defined as ob- 
jects with magnitude 18.0 < r < 19.5 that are identified 
as pointlike by the SDSS photometric pipeline |5ll). and 
a dust reddening map of E{B — V) from Ref. |8J]. These 
maps can be substituted in place of the LRG map Sug/fig 
as XiflG in Eq. (PS|l and the cross-spectrum and "bias" 
determined. The biases obtained using these contami- 
nant maps are shown in Table IVTI 

A crude estimate of how this contamination trans- 
lates into contamination of the galaxy-convergence power 
spectrum Cf^ is provided by performing an unweighted 
least-squares fit of the LRG density to the reddening and 
stellar maps, 



Ug rii. 



ci + residual, (49) 



over the 296,872 SDSS pixels. The fit coefficients are 
CE = 0.62, c^ = -0.0092, and ci = -0.017. We have 
shown in Table IVTI the spurious contribution Abg to the 
bias resulting from stellar and reddening contamination 
if one assumes these fit coefficients. 



2. Microwave foreground template test 

Since the lensing map v is a quadratic function of tem- 
perature, and the Galactic foregrounds are not corre- 
lated with the primary CMB, Eq. (|40|l is applicable to 
Galactic foregrounds. The contamination A(v"'') can be 
obtained straightforwardly since the right-hand side of 
Eq. (|40|l can be evaluated by substituting in the fore- 
ground maps as AT". The difficult step is to construct 
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a good foreground map AT"; here we use external tem- 
plates to avoid any possiblity of spurious correlations of 
the templates with the WMAP data (either CMB signal 
or noise). 

The Galactic foregrounds that must be considered in 
producing a template at higher (W band) frequencies are 
free-free and thermal dust emission; at lower frequen- 
cies (Q and V bands) an additional component is present 
whose physical origin remains uncertain but which may 
include hard synchrotron emission J83l | or spinning or 
magnetic dust 99, 100]. We have used Model 8 of Refs. 
[84, 85] for thermal dust, and the Hcu line radiation tem- 
plate of Ref. [Sy] re-scaled using the conversions of Ref. 
J83J for free-free radiation. There are no all-sky syn- 
chrotron templates at the frequencies and angular scales 
of interes t ( the Haslam radio continuum maps at 408 
MHz 1^ 1^3 ) frequently used as a synchrotron template 
for CMB foreground analyses, have a 50 arcmin FWHM 
beam and hence do not resolve the I ~ 200-400 scales 
used for lensing of the CMB). Nevertheless, inclusion of 
the low-frequency component (w hatev er its origin) is not 
optional, and so we follow Ref. |lO0l ] in modeling it as 
proportional to the thermal dust prediction of Refjs3 
multiplied by T|„^j using the coefficients of Ref. jlOOJ . 

As a test for contamination, we have substituted these 
foreground templates for the true CMB maps, run them 
through the lensing pipeline, and derived bg estimates 
by correlating against the true LRG map; the results 
are shown in the "foreground" column of Table |n] The 
typical contamination due to foregrounds is clearly very 
small (bias error of a few times 10~^), and thus is neg- 
ligible even if the foreground amplitude has been un- 
derestimated by an order of magnitude (the error on 
bg scales as the foreground amplitude squared). This 
is not surprising: in the relatively clean regions of sky 
used for this analysis, the galactic foreground temper- 
ature anisotropy is roughly 2 orders of magnitude (in 
amplitude) below the CMB temperature at the degree 
angular scales. Therefore a quadratic statistic, such as 
the lensing estimator, should be ^ 4 orders of magni- 
tude smaller than the foregrounds (again, in amplitude). 
Thus when v (foreground) is correlated against the galaxy 
map, the correlation that one expects from chance align- 
ments of foregrounds and galaxies is roughly 4 orders of 
magnitude less than v(CMB) (although a much greater 
correlation could exist if the galaxy map were also con- 
taminated by foregrounds, e.g. dust extinction). This is 
in contrast to the point sources, which are highly local- 
ized objects that become more and more dominant when 
we consider higher-order statistics such as the lensing es- 
timator V. 



VII. DISCUSSION 

In this paper, we have carried out an initial search 
for weak lensing of the CMB by performing a lensing 
reconstruction from the WMAP data and correlating the 



resulting lensing field map with the SDSS LRG map. We 
do not have a detection, however our result bg = 1.81 ± 
1.92 (Ict) is consistent with the bia s bn ~ 1.8 obtained 
from the LRG clustering autopower [601 ]. 

The main purpose of this analysis was to identify any 
systematics that contaminate the galaxy-convergence 
correlation at the level of the current CMB data. The 
good news is that our result for the bias is reason- 
able, suggesting that such systematics are at most of 
the order of the statistical errors. We have also found 
that the Galactic foregrounds are a negligible contam- 
inant to the lensing signal (again, at the level of the 
present data). The bad news primarily concerns extra- 
galactic foregrounds: a significant amount of solid an- 
gle - 21% of Kp05nS10 - was lost due to point source 
cuts that are necessary to avoid spurious power (at least 
in Q band), and we have no assurance that significant 
point source or tSZ contamination of the lensing sig- 
nal does not lie just below the threshold of detectabil- 
ity. The extragalactic foreground analyses of Sec. IVI D II 
based on the frequency dependence of the signal and 
the galaxy-temperature correlations yielded only a weak 
constraint on the synchrotron point source contamina- 
tion, Abi^'^\TT) = 0.73 ± 1.18, and essentially no use- 
ful constraint can be derived for tSZ using the first-year 
WMAP data. Our constraints based on Cf for the point 



sources are more stringent, Ab 



(PS) 



-0.14 ±0.51, but 



these assume Poissonianity of the sources, which must 
break down at some level. The point source and tSZ 
issues will become even more important as future ex- 
periments probe lensing of the CMB using higher-/ pri- 
mary modes, where point source and tSZ anisotropics 
contribute a greater fraction of the total power in the 
CMB, and precision cosmology with lensing of the CMB 
will require a means of constraining these contaminants 
in order to produce reliable results. Because in the real 
universe the extragalactic foregrounds will not be exactly 
Poisson-distributed, the frequency (in)dependence of the 
lensing signal will in principle provide the most robust 
constraints on the contamination. In this paper, we were 
unable to obtain useful constraints this way because of 
the limited range of frequencies on WMAP (all in the 
Rayleigh- Jeans regime where the tSZ signal has the same 
frequency dependence as CMB) and the low statistical 
signal-to-noise. Both of these problems should be alle- 
viated with high-resolution data sets covering man y fre - 
quencies, e.g. as expected from the P7aiicJf satellite 'ill]. 

The large solid angle that was lost to point source cuts 
in the analysis presented here resulted from the need 
to remove "artifacts" in the v map that occur around 
point sources. One approach to this problem would be 
to try to devise a CMB lensing reconstruction technique 
that works with complicated masks. Alternatively, one 
could compute the galaxy density-CMB temperature- 
CMB temperature bispectrum Bf ^ j , rather than try- 
ing to use the lensing field (a quadratic function of CMB 
temperature) as an intermediate step; this way one could 
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mask out only the point source itself and not a 2 degree 
exclusion radius around it. The bispectrum approach 
carries the additional advantage of retaining angular in- 
formation about the foregrounds; this information may 
be useful for separating lensing from the kinetic SZ and 
patchy-reionization anisotropies that have no frequency 
depedence but can still contaminate lensing if small- 
scale information is used 79, 94, 101]. The bispectum 
Bi I I may therefore be of particular interest for lens- 
ing analyses of high-Z ex perim ents such as the Atacama 
Pathfinder Exp eriment |112| . the Atacama C osm ology 
Telescope [ulf, and the South Pole Telescope |lT^ . 

In summary, this paper represents a first analysis of 
lensing of the CMB using real data, and should not be 
regarded as the last word on the methodology. We have 
identified extragalactic foregrounds (point sources and 
tSZ) as the most worrying contaminant to the lensing sig- 
nal in the WMAP data; the point sources, if unmasked, 
dominate the power spectrum of the reconstructed con- 
vergence if the Q band data are used, but this effect is 
suppressed at higher frequencies. We have shown that in 
the current data, the Galactic foreground contribution 
is negligible, and the contamination from point sources 
and tSZ in the galaxy-convergence cross-spectrum is at 
most of order the signal (although we have no detection 
of contamination). Like the search for the CMB lens- 
ing signal (and its eventual use in precision cosmology), 
stronger statements about the foreground contamination 
must await higher signal-to-noise data at a wide range of 
frequencies. 
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APPENDIX A: NON-ISOLATITUDE SPHERICAL 
HARMONIC TRANSFORM 

The non-isolatitude spherical harmonic transform 
(SHT) is used in our cross-correlation analysis. The SHT 
operations on the unit sphere transform between real- 
and harmonic-space representations of a function: 

r(n,) = J2'iZo Y!m=~i Ti„,Yi^{-h.i), (synthesis) 
Si^ = E^V' YCm{■^^)S{i^i)■ (analysis) 

(Al) 
For high-resolution data sets, this operation is usually 
performed using an isolatitude pixelization, i.e. one in 
which the pixels are positioned on curves of constant co- 
latitude 9. This situation allows the colatitude {9) and 
longitude (cf) parts of the spherical transform to be per- 
formed independently, resulting in an overall operation 
count scaling as 0{Np ) |l02lll03l ]. While this approach 
works, and has contributed remarkably to the popular- 
ity of isolatitude pixelizations such as HEALPix, there 
are reasons to maintain the flexibility to use any pixels. 
For example, in simulations of gravitational lensing of the 
CMB, we need to produce a simulated lensed map, and in 
general a set of pixels that are isolatitude in "observed" 
coordinates (e.g. HEALPix) maps onto a non-isolatitude 
grid on the primary CMB. We note that for other analy- 
ses there may be other reasons to consider more general 
pixelizations, whi ch preserve desired pro pert ies such as 
conformality |l04j or maximal symmetry |l05l ] . This Ap- 
pendix describes our non-isolatitude SHT algorithm. 



1. The method 

We consider first the SHT synthesis. Our first step is to 
perform a latitude transform using associated Legendre 
polynomials on a set of points equally spaced in 9 (the 
"coarse grid"): 



Tm( — Vr) = 2_^ TimYlmi 
/— |7n| 



L 



7r,0 = O). 



(A2) 



This procedure is performed for integers in the range 
< a < L, —/max < "T' < /max- [Here L is an inte- 
ger satisfying L > /max, which we require to be a power 
of two times a small odd integer. The first requirement 
ensures that Eq. IJA2|I over-Nyquist samples the varia- 
tions in Tm {9) , the second ensures that the FFT is a fast 
operation.] It requires a total of 0{l'^^^L) operations and 
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is the most computationally demanding step in the trans- 
form. The spherical harmonics are computed as needed 
using an ascending recursion relation. The standard re- 
cursion relation is used to generate the associated Legen- 
dre functions; we speed up the transform by a factor of 
two over the brute-force approach by taking advantage of 
the symmetry/antisymmetry of the spherical harmonics 
across the equatorial plane. 

The next step is to refine the coarse grid, which has a 
spacing of tt/L in 9, to a "fine grid" with spacing tt/L' 
where L' > L. We do this by taking advantage of the 
band-limited nature of the spherical harmonics. Any lin- 
ear combination of spherical harmonics of order I < /max 
can be written as a band-limited function: 



T,„(0) = 



l — \7n\ 



^ Ini ^ Ini \ 



9,0) = 



E 



Or, 



ijiO 



(A3) 



We may determine the coefficients Cm.n via a fast Fourier 
transform (FFT) of length 2L, so long as the left-hand 
side has been evaluated at the points 9 = na/ L for inte- 
gers —L < a < L. (We use parity rules to compute the 
left hand side at negative values of 6*.) By applying an 
FFT of length 2L' to the C,n,n, we then obtain T„,{9) at 
values of = na/L' for integers < a < L' . What has 
been gained here is that we have performed the associ- 
ated Legendre transform on 2L' points, but the expensive 
evaluation of the associated Legendre polynomials has 
only been required at 2L points. If L and L' are powers 
of two, then the FFT process requires 0{lmaxL'\ogL') 
operations. 

The third step is an FFT in the longitude direction 
to obtain T{9 = ira/L' , (p = n-f/L'), where a and 7 are 
integers. This process has become standard in isolatitude 
SHT algorithms; in its full glory, it is given by: 



iratj) 



T[a,7]^T(0=-^7r,<^=^7r)= £ r„,(0)e 



(A4) 
At the end of this step, we know the real-space value of 
our function on a fine equicylindrical projection (ECP) 
grid of spacing tt/L' . The FFT operation in this step 
requires 0{L'^\ogL') operations. The total operation 
count of transforming onto the ECP grid is 0(/max) ^^'^ 
is dominated by the associated Legendre transform, so 
long as L' /l^sx < \//max/log/max- What is important to 
note is that, particularly if /,„ax is large (for a future CMB 
polarization experiment we have to consider multipoles 
up to roughly /max ^ 4000), we can sample our function 
in real space at several times the Nyquist frequency at 
no additional computational cost. This is exactly what 
is required in order to successfully interpolate the values 

The final step is the interpolation step. For each point 
hj, we identify the coordinates in the ECP grid: (we 
suppress the j index here for clarity) 



a + 5a — L' — ; 7 -f (5^ 



L'- 



(A5) 



where a and 7 are integers and the fractional parts < 
(5q,,(5^ < 1. a 4iir^-point, two-dimensional polynomial 
interpolation is then computed: 



Wy{5^)T[a+^,"f+v], (A6) 



where the weights Wp{5) are computed by Lagrange's for- 
mula: 
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Wp((5) = 



(-l)^-P 



[K-p)\[K-\ + p)\{6-p) 



n (^ 



(A7) 

The weights for both the a and 7 directions may be 
evaluated in a total of 0(IC) multiplications and divi- 
sions if the factorials have been pre-computed, so that 
for high-order interpolations the dominant contribution 
to the computation time in interpolation comes from the 
multiplications in Eq. (|A6|) rather than from computa- 
tion of the weights. 

Note that the "analysis" operation of Eq. (|A1() is the 
matrix transpose operation of the "synthesis" (if we view 
the pixelized T{\\.i) and the harmonic-space Tim as vec- 
tors), and that all of the steps outlined above are linear 
operations. Since the transpose of a composition of op- 
erations is the composition of the transposes in reverse 
order, we can simply use the transposes of these opera- 
tions in reverse order to compute an SHT analysis using 
the same number of operations. 

The above SHT algorithm generalizes easily to vector 
spherical harmonics, for which we use the basis: 

1 






1 



vri,„(n). 



fix Vy;„i(n). 



(A8) 



2. Interpolation accuracy 

The order K and step size 7r/i' of the interpolating 
polynomial is determined by a balance of computation 
time, memory usage, and accuracy. The computational 
cost of the interpolation is 0{K'^Npix) and the memory 
usage is 0(L' ). Hence it is important to understand how 
interpolation accuracy is related to K and L' . 

The error in polynomial interpolation of T{9, 0) can 
be determined from the error in interpolation of a band- 
limited function. We note that the error in interpolation 
of a given Fourier mode f{9, 
by: 

finterp{0, <P) 



C™^„e*("^+™*'-' is given 



I{e.^) 



/ Tfi \ 1 r / ix \ 



where the v function is: 

K 



(A9) 
(AlO) 
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Error in polynomial interpolation 




1e-07 n 



1e-08 



0.4 0.6 

Phase advance, v/jt 

FIG. 13: The maximum fractional error jwmax| from Eq. l|A12fl 
as a function of the interpolation order K and the sampling 
rate. From top to bottom, the curves correspond to 7^ = 1 
to K = 10. Note that tp = it for data sampled at the Nyquist 
frequency (i.e. L' = Zmax), whereas tp = n/2 for data sampled 
at twice the Nyquist frequency (L' — 2Zinax), etc. Accuracy 
can be improved via either fine sampling (small i/>) or high- 
order interpolation (high K). 



(Note that ip represents the phase advance of our Fourier 
mode per grid ceU.) The v function determines the frac- 
tional error in a given Fourier mode. It is most easily 
evaluated by noting that (since polynomial interpolation 
is exact for a constant) we have v(5, 0) — 0, and deriva- 
tive: 



dv 



dip 



n 



K 



{2K-1) 



II 



-i-iji\2K-l 



(All) 



[Here we have used Eq. IJA7|I and applied the binomial 
theorem.] We note that the product is maximized at 
5 = 1/2 and that using trigonometric identities the ex- 
ponential term can be simplified. Since v{5, 0) = 0, we 
have that |t;| cannot exceed the integral of \dv/d'4>\ from 
to V: 



\v{5,^)\< 



{2K)\ 



A^K\{K ~ 1)! 



'sin^^-i^d^'. 
2 



(A12) 



This bound is plotted in Fig. 1131 



APPENDIX B: PRECONDITIONER 

One of the steps in the cross-spectrum estimator re- 
quires the solution of the linear system: 



CpY 



(Bl) 



for y. Here Cp is the prior covariance matrix for the 
LRGs and is given by Cp = Sp -I- N , where Sp is the signal 



prior and N is the noise. Since the matrix Cp is too large 
(iV(iV+l)/2 elements where the dimension N = 296, 872) 
to store in memory, Eq. (|B1|I must be solved using itera- 
tive methods. We have used a preconditioned conjugate- 
gradient (PCG) method. The PCG method 10&1l07J 
with preconditioner E produces successive estimates y^*^ 
for the solution to Eq. (|B1|I using the equation: 



rii) 



,(i-l) 



j.(*-l)T£j.(-*-l) 



pWCpP 



,(■0 



(0 



(B2) 



where the residuals are defined by r^*' = x — Cpy'-*' , 
the search directions p^'^ are chosen according to 



and 



j(») ^ Er'^"^) 



j.(j-l)T£j.(j-l) 

r(*-2)TEr(«-2)' 



.(«-!) 



(B3) 



The initial conditions are y'"' = and p'"'^' = Ex. 

The choice of preconditioner strongly affects the rate 
of convergence of the PCG algorithm; ideally we have 
a preconditioner for which Eu can be rapidly computed 
for given u, and for which E w C""'^. I n usi ng PCG for 
CMB power spectrum estimation, Ref. |l06j used a pre- 
conditioner based on (approximate) azimuthal symmetry 
of the Galactic plane cut. Unfortunately, the SDSS sur- 
vey region lacks any such symmetry. Therefore we use a 
preconditioner that works independently of any symme- 
tries of the survey region. The strategy is to use a two- 
scale preconditioner: the low-resolution scale (/ < IspUt) 
is solved by brute-force matrix inversion in harmonic 
space, while the high- resolution scale {I > IspUt) is essen- 
tially un-preconditioned and the burden of convergence 
lies on the conjugate-gradient algorithm. This is useful 
because the large condition ratio of a typical prior ma- 
trix Cp comes almost entirely from a few large eigenvalues 
corresponding to the low / modes. The condition ratio is 
dramatically improved by suppressing these eigenvalues 
with the two-scale preconditioner. 

The preconditioner is obtained by cutting the power 



spectrum at IspUt and defining the I 
M: 



split 



X I 



split 



matrix 



\/{Sp,l - ■S'p,i,pHt)('S'p,i' - ■S'p.i.pHt) 

>^T.^jyL{nj)Y,„Anj). (B4) 



Here I and I' are in the range from to Igput — 1, j is 
a pixel index, and flj is the area of the jth pixel. Note 
that a red prior spectrum (i.e. Spj > Sp^i^^^.^) has been 
assumed, as appropriate for the angular power spectra of 
galaxies. We next define the G matrix: 



Gt, 



hn.l' 



\J{Sp.l - Sp^i^^,^^){Sp,y - Sp^i^^,^^) 



5. 



p.hplit 



x\M-^' 



Im.l'' 



(B5) 
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Then the preconditioner is: 



Ein — 



i L-iO-ij 



u u^ u u 'J 






(B6) 



Using the Woodbury matrix inverse formula, the ma- 
trix E of Eq. IJB6|I can be recognized as the exact inverse 
of the matrix: 






hplit — ^ 

E 

1=0 



{Sp.i 



^P-lsplit) 



m^ — l 



Yir,.{n,)Y:^{h,). 



(B7) 



That is, E is the exact inverse of a convolution that has 
power spectrum Sp^i for I < Isput and white noise with 
power Sp,i^ ,.j for I > Igput- Thus E is a good approxi- 
mation to Cp if Ispiit is chosen to be the multipole where 
the signal and noise have similar power. In this case, 



the convergence of the PCG iteration is very rapid. In 
practice (see below) we choose Isput to be somewhat less 
than this in order to speed up multiplication of a vector 
by E, and accept that more iterations will be required for 
convergence. 

Note that L need only be computed once for each sky 
cut and prior power spectrum; it can then be stored 
and used for many Eu operations. The Eu operation 
is then reduced to a spherical harmonic transform of cost 
0[l%ut) + 0{K^Npix), a matrix-vector multiplication of 
computational cost 0{l'^ i^^), and another spherical har- 
monic transform. Since we use the preconditioner once 
for every Cp operation, we gain speed by increasing Igput 
until the most expensive part of the E operation is of com- 
parable cost to the Cp operation (O(^max))- This suggests 
that we set I spin ~ ^max- There are however practical 
limitations on Ispia '■ the size of the matrix in memory is 
0{ltpiit)i and the computational cost of obtaining G is 
0{l^ m). Thus the best value of IgpUt is generally some- 



what less than Zmax; we have used I 



split 
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